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INTRODUCTION 


The inclusion of body forces has received a good deal of attention in boundary element research. 
The consideration of such forces is essential in the design of high performance components such as 
fan and turbine disks in a gas turbine engine. Due to their critical performance requirements, optimal 
shapes are often desired for these components. The boundary element method (BEM) offers the 
possibility of being an efficient method for such iterative analysis as shape optimization. 

A survey of efforts in the area of sensitivity analysis in BEM was given by Mota, Soares and Choi 
[1]. The shape sensitivity using a finite-difference formulation was given by Wu [2] and using the 
implicit-differentiation formulation by Barone and Yang [3], Saigal et al. [4-6], and Rice and 
Mukheijee.* Mukherjee and Chandra [7] presented a BEM sensitivity formulation for materially 
nonlinear problems. The treatment of body forces for sensitivity analysis has not received much 
attention. 

In this paper, the implicit-differentiation of the boundary integral equations [8] is performed to 
obtain the sensitivity equations. The body forces are accounted for by either the particular integrals 
[9,10] for uniform body forces or by a surface integration [1 1] for non-uniform body forces. The 
corresponding sensitivity equations for both these cases are presented. The validity of present 
formulations is established through a close agreement with exact analytical results. 


*Rice, J.S. and Mukherjee, S., "Design Sensitivity Coefficients for Axisymmetric 
Elasticity Problems by Boundary Element Methods", private communication. 
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BOUNDARY ELEMENT ANALYSIS EQUATIONS 


Including the effect of temperature variation <(), the stress tensor Ojj is given in equation (1), and the 

equation of equilibrium is given in equation (2). Starting with a weak statement of equation (2) and 
using the divergence theorem twice the integral equation (3) is obtained. Assuming steady state 
condition, using the divergence theorem, and applying Green's second function leads to equation (4) 
where the thermal effects have been reduced to a boundary integral form. 

F vF F 

Ojj — 6;; + —————— T^iiad) (1) 

,J (1+v) J (l+v)(l-2v) (l-2v) ,ja * 


a ijj + F i 


= 0 


Uj(p) = - J r T ij (p,q)u j (q)dr + J Ug(p,q)tj(q)dr + J^U ij (p,q)F i (q)dQ - 

TTvrf UijCp.q^dn 

(1— 2v) Jq 

Uj(p) = -J^T i j(p,q)u j (q)dr + |^U ij (p,q)t j (q)dr + 


«ii+vLf{±r n r .^_ r m 

87i(l-v)Jrl RV dnj dnj 


dr+ J^U ij (p,q)F i (q)dQ 


(2) 


( 3 ) 


( 4 ) 


E, v, and a are the modulus of elasticity, Poisson's ratio, and coefficient of thermal expansion, 
respectively. T- and Uy are the fundamental (Kelvin) solutions for traction and displacement, 
respectively, p and q are the load point and the field point, respectively, and R is the distance 
between these points. u v t; , and F ; are the components of displacement, traction and body force, 
respectively. 

GRAVITATIONAL AND CENTRIFUGAL FORCE SENSITIVITY 

If the total displacement is written as a sum of a complementary and a particular integral component 
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as in equation (5), then in the absence of temperature variation, the last two terms in equation (4) 

drop out giving equation (6). Discretizing the boundary using boundary elements with 

displacements and tractions interpolated as shown in equation (7), we get the matrix relationship 

(8). Substituting from equation (5) into equation (8), we obtained a relationship in equation (9) 

including the effect of particular solutions due to body forces. Particular solutions {uP} and {tP} 

were given by Banerjee and co-workers [9- 10] > at SUNY - Buffalo. Implicit differentiation of equation (9) 

with respect to the design variable Xj^ results in the sensitivity equation (10). The contribution of 


the body forces is included in the vector {fP} given in equation (11). 

Uj = u® + uf (5) 

Ui(p) =-| T ij (p,q)u J c (q)dr + f Ujj(p,q)t?(q)dr (6) 

Discretizing equation (6) using interpolation functions for displacements and tractions 
u = [H]{u} ; t = [H]{t} (7) 

The matrix form of equation (6) is obtained as 

[F]{u c } = [G]{t c } (8) 

Substituting equation (5) in equation (8) 

[F]{u> = [G] {t} + [F] {u p } - [G] {t p } (9) 

Differentiating with respect to the design variable, X L 

[F]{u}, l = [G], L {t} + [G]{t}, L -[F], L {u} + {?} (10) 


Where 


[?} = [FKu p }, L + [FJ, L {u p } - [G]{t p }, L - [G], L {t p > 


( 11 ) 



The superscripts c and p refer to the complementary and particular soultions, respectively. [ H ] is 
a matrix of interpolation functions. ( ), L denotes the derivative of ( ) with respect to the design 
variable X L . 


THERMOELASTIC SENSITIVITY 


For the case of temperature variation $ and temperature gradient <j), n , the term with volume integral 

in equation (4) drops out. Then using the interpolation given in equation (7), we get the matrix 
relationship given in the equation (12). Implicit -differentiation of equation (12) leads to an equation 
similar to equation (10) but with a different definition for vector {fP}. This relationship is given in 
equations(13) and (14). The matrix [V] involves thermoelastic kernels which include elliptic 
integrals of the first and the second kind. The present sensitivity analysis requires derivatives of 
these elliptic integrals which are easily determined through chain rule of differentiation. 


[F] {u} = [GKO + [V]{T} 


( 12 ) 


[F]{u}, L = [G}, l {0 + [G](0 , l - [FLiiiO + {f} 


(13) 


{( p } = [V], L {T} + [V]{T}, L (14) 

(T) is the vector of nodal temperatures. 
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SEMI- ANALYTICAL SENSITIVITY FORMULATION 


The sensitivities can now be obtained using equations (10) and (11) for centrifugal and gravitational 
body forces, and using equations (13) and (14) for thermal body forces. We, however, need to 
determine sensitivity matrices such as [F},^ and [G],jj and sensitivity vectors such as {uP},l and 

{tP},L. In the semi- analytical approach, the design variable Xl is first perturbed by an amount 

AXj^. The system matrices [F(Xl+AXjJ] , [GCXl+AX^)], etc., are generated based on the new 

geometry. The sensitivities are then simply obtained using forward-difference relationships shown 
in equations (15) and (16). It is noted that the sensitivity results will depend on the perturbation step 

size AXj^. However, this step will result in substantial simplification of the implementation of the 
sensitivity algorithm. 


[F(X l +AX l )] -[F(X l )] . 

[f] ' l = ES ’ etc 


( 15 ) 


{u p }, L = 


[u p (X l +AX l )1 - lu p (X L )l 

ax l 


etc 


( 16 ) 


FULL ANALYTICAL SENSITIVITY FORMULATION 

For the full-analytical approach, the sensitivity matrices and vectors are directly calculated from then- 
analytical expressions given in equations (17) and (18). These expressions, however, need the 
sensitivities of geometry quantities such as x.^, y,L » n >L> etc - Th e i n itia.I geometry is first used for 
solution of vectors (u) and (t) in equations (8) or (12). This geometry is then changed through a 
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perturbation AXj^ of the design variable. Only the geometry sensitivities are then calculated using 

forward-difference approximation. These geometry sensitivities are needed for evaluating terms in 
equations (17) and (18). 


N 1 

[F], l = I J„ { [t*],5 H] J + [t*] T [ H] J, L } 



The superscripts T and * refer to the transpose of the matrix and the fundamental solutions, 
respectively. 

3-D Centrifugal Loading Particular Integral Sensitivities: 

ul L =( 2c 1 (xx tL +yy iL )+2c 2 zz tL )x+(c 1 (xx+yy)+c 2 zz)x tL 
u 2, l =( 2c 1 (xx >L +yy iL )+2c 2 zz >L )y+(c 1 (xx+yy)+c 2 zz)y ;L 
u3, L =2c 3 (xx iL +yy >L )z+c 3 (xx+yy)z (L 

C!= -pto 2 (5A,+4ji)/4(X+u) 
c 2 = -pco 2 (i/8(X+2(i)/(X+ji) 

c 3 = pco 2 /8/(X,+2|i) (18) 


X, fi are the Lame's constants. x L denotes the derivative of coordinate x with respect to the design 
variable X L and y L and z L have similar definitions.The corresponding traction sensitivities can be 
found using the constitutive relation. 
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SINGULAR TERMS IN SYSTEM SENSITIVITY MATRIX 


For determining terms in [F],l an extension of rigid body technique used for singular terms in [F] 
is used. This extension is based on the fact that the sensitivities corresponding to rigid body 
displacements and tractions are zero leading to a row-sum type property for [F],^. Thus from 

equation (10), the singular terms for 3-D can be obtained as given in equation (21). For 2-D, these 
terms are similarly obtained. For axisymmetric case , a rigid body motion in Z - direction and an 
inflation mode in the radial direction are used. 


"f u rigid^xi _ 


'1 1 


o' 

0 


i 

1 


0 

0 

> * ^ u rigid^X2"”' < 

1 > 

1 


0 

.0 . 


. 1 . 


(19) 


{t}rigid — {0} 


( 20 ) 


n-1 n 

X Fjj, L = 0 ; Fjj, L = 0 ; for each i 

j= 1,3,*** j = 2,4,- 


( 21 ) 
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STRESS SENSITIVITY RECOVERY 


The solution of sensitivity equations yields the boundary traction sensitivities only. The stress 
sensitivities at other locations can be obtained directly using differentiated elasticity equations 
describing stress-Strain relationships. These relations are given in equation (22) for the 
axisymmetric case. 


<*ir' 


r 2 2 0 - 

n r ti z 2n r n z 


^ <^11 ' 

^zz 

> = 

2 2 0 

n z i\ 2n r n z 

< 

a 22 y 

^rz j 


_ n r n z -n r n z (i\ 2 -n z ) _ 


- a 12- 


r^, n z are the components of the outward normal in the r and z directions, respectively. 


C 22T“ 77"“T a ll,L + ^r( e 22,L+ Ve 8G,d 

U-V) d_ v 2 ) 

_ J, L 1 

e 22,L ~~~2 U2 £ + J U2 >$ L 


eee,L = ~r(u r>L r - u r r L ) 


( 22 ) 


U 2,^L d enotes the mixed derivative with respect to the dimensionless coordinates % and the design 
variable X L . 
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NUMERICAL EXAMPLES 


The above formulations were applied to a series of selected examples to determine the design 
sensitivities for displacements and tractions. These examples include: (a) a rotating circular disk of 
constant thickness with a central hole shown in Figure 1 analyzed using two-dimensional elements, 

(b) a rotating circular disk with hyperbolic varying thickness and with a central hole shown in 
Figure 2 analyzed using axisymmetric elements, (c) a hollow cylinder under plane strain shown in 
Figure 3 subjected to pressure and temperature change and analyzed using axisymmetric elements, 

(d) a pressurized hollow cylinder under temperature variation shown in Figure 4 and analyzed using 
axisymmetric elements, (e) a solid circular bar (.Figure 5) under self-weight analyzed using three-dimensional 
elements, and (f) a rotating circular disk (Figure 6) analyzed using three-dimensional elements. For examples 

without temperature variation the material data used were E=30xl0^ psi, v=0.3; and for examples 

with temperature variation the material data used were E=1 psi, v=0.3, and a=0.02 /°F. The results 
obtained from the present formulations were compared with the exact solutions to check these 
formulations. For exact sensitivities the elasticity solutions were first expressed in terms of the 
design variable and then differentiated with respect to this design variable. A good comparison of 
the present results was seen from the results presented in the following pages. 
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Design Sensitivity Analysis of a Rotating Circular Disk. 


Sensitivity 


Radial Displacement 

x (10 ' 3 ) 

Traction (x 10 3 ) 


Location* 

Analytical 

This 

Analytical 

This 



Study 


Study 

A 

6.7678 

6.7677 

4.8061 

4.7625 

B 

6.5929 

6.5929 

8.7301 

8.6993 

C 

2.8329 

2.8330 

3.6905 

3.6923 

D 

2.6684 

2.6685 

3.8066 

3.8103 

E 

6.7678 

6.7678 



F 

6.7678 

6.7676 




Location coordinates (r, z) in inches are: A(4.333, 0.), B(5.061, 0.), C(17.061, 0.), D(19.0, 0.), 


E(3.864, 1.035), F(3.967, 0.522) 



”*“4IN 


16 IN 


FIGURE 1 
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Design Sensitivity Analysis of a Rotating Disk 
with Hyperbolic Varying Thickness 




SENSITIVITY 



RADIUS 

(inch) 

DISPLACEMENT X 10 3 

HOOP STR. X lO- 3 


Exact 

Mesh A 

Mesh B 

Exact 

Mesh B 

4.00 

2.05931 

1.9987 

1.9986 

5.07485 

3.9470 

4.25 

2.05819 

- 

1.9887 

5.66207 

5.4810 

4.50 

2.04760 

1.9725 

1.9725 

5.99760 

5.6433 

4.75 

2.02989 

- 

1.9522 

6.16066 

5.8450 

5.00 

2.00685 

1.92% 

1.9294 

6.20577 

5.8609 

5.25 

1.97984 

- 

1.9046 

6.17095 

5.8355 

5.50 

1.94995 

1.8789 

1.8785 

6.08307 

5.7574 

5.75 

1.91803 

- 

1.8512 

5.%129 

5.6542 

6.00 

1.88478 

1.8236 

1.8229 

5.81938 

5.5318 

6.25 

1.85078 

- 

1.7940 

5.66731 

5.4024 

6.50 

1.81653 

1.7653 

1.7645 

5.51231 

5.2668 

6.75 

1.78246 

- 

1.7348 

5.35968 

5.1372 

7.00 

1.74893 

1.7060 

1.7052 

5.21325 

5.0076 

7.25 

1.71628 

- 

1.6759 

5.07583 

4.8915 

7.50 

1.68481 

1.6482 

1.6474 

4.94945 

4.7787 

7.75 

1.65480 

- 

1.6198 

4.83557 

4.6831 

8.00 

1.62650 

1.5945 

1.5937 

4.73521 

4.5934 

8.25 

1.60016 

- 

1.5692 

4.64908 

4.5222 

8.50 

1.57602 

1.5474 

1.5466 

4.57767 

4.4589 

3.75 

1.55428 

— 

1.5261 

4.52126 

4.4141 

9.00 

1.53517 

1.5089 

1.5081 

4.48003 

4.3782 

9.25 

1.51889 

- 

1.4928 

4.45403 

4.35% 

9.50 

1.50565 

1.4811 

1.4803 

4.44327 

4.3511 

9.75 

1.49564 

- 

1.4710 

4.44769 

4.3540 

10.00 

1.48906 

1.4660 

1.4652 

4.46718 

4.4120 


Note: Exact solution is for the assumption of plane stress. 
Mesh A: 15 element model; Mesh B: 30 element model 


I 
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DESIGN VARIABLE: INNER RADIUS, R. 

I 


FIGURE 2 



Design Sensitivity Analysis of a Plane Strain Hollow 
Cylinder Under Pressure and Temperature Variation 


SENSITIVITY 


RADIUS 

(inch) 

RADIAL DISPLACEMENT 

RADIAL STRESS 

CIRCUMFERENTIAL 

STRESS 

Exact 

This 

Study 

Exact 


Exact 


3.0 

3.6305 

3.6315 


0.000 

wm 

0.5956 

3.5 

3.6694 

3.6712 


-0.1074 

■ss 

0.6993 

4.0 

3.6168 

3.6174 


-0.1186 

BH 

■ 

4.5 

3.5302 

3.5309 


-0.0975 

wBm 

1 £ 

5.0 

3.4341 

3.4350 


-0.0654 


0.6594 

5.5 

3.3396 

3.3402 


-0.0318 

0.6264 

0.6265 

6.0 

3.2513 

3.2523 

0.000 

0.000 

0.5955 

0.5956 



















Pressurized Hollow Sphere Under Temperature Variation 


RADIUS 
(inch) 


RADIAL DISPLACEMENT 


RADIAL STRESS 


CIRCUMFERENTIAL 

STRESS 



Exact 

HBH 

Exact 

HSE9HI 

Exact 


1.00 

1.5126 

1.4992 

0.000 

0.0003 

HI 

1.495 

1.125 

1.758 

1.7448 

-0.754 

-0.618 


1.886 

1.25 

1.8563 

1.8436 

-0.901 

-0.927 

111311 

1.910 

1.375 

1.8985 

1.887 

-0.815 

-0.791 

HHi 

■ESS 

1.50 

1.9247 

1.9144 

-0.653 

-0.642 

1.801 


1.625 

1.9526 

1.9433 

-0.471 

-0.466 

1.715 

■HI 

1.75 

1.9896 

1.9813 

-0.297 

-0.290 

1.634 

1.630 

1.875 

2.0384 

2.0311 

-0.139 

-0.128 

1.562 

1.561 

2.00 

2.0996 

2.0952 

0.000 

0.0003 

1.500 

1.496 


















Solid Circular Bar under Self-weight 


Location Sensitivity of displacements in Z direction (10^) Sensitivity of Stress in Z direction 


X =3.536 
y =3.536 

Exact 

Full- 

Analytical 

Semi- 

Analytical 

Exact 

Full- 

Analytical 

Semi- 

Analytical 

A(z=0.0) 

0.0 

0.0 

0.0 

2318.4 

2320.4 

2351.6 

B(z=5.0) 

6.7620 

6.7613 

6.7614 

1738.8 

1747.6 

1758.3 

C(z=10.0) 

11.592 

11.589 

11.589 

1159.2 

1188.7 

1189.8 

D(z=15.0) 

14.490 

14.488 

14.488 

579.6 

587.6 

589.4 

E(z=20.0) 

15.456 

15.448 

15.448 

0.0 

0.0 

0.0 



4 


20 IN 




FIGURE 5 



Three Dimensional Rotating Circular Disk 


Location Displacement Sensitivity 11^(1 O’ 3 ) Radial Stress Sensitivity (10 4 ) 

y=0.0 

z=1.5 



Exact 

Full- 

Analytical 

Semi- 

Analytical 

Exact 

Full- 

Analytical 

Semi- 

Analytical 

X 

II 

O 

0.6768 

0.6775 

0.6735 

0.0840 

0.1907 

0.2023 

O 

iri 

II 

X 

0.6613 

0.6610 

0.6571 

0.8564 

0.8413 

0.8511 

X 

II 

o\ 

o 

0.6216 

0.6212 

0.6179 

0.9650 

0.9674 

0.9644 

x= 7.0 

0.5768 

0.5765 

0.5737 

0.9018 

0.8966 

0.8945 

X 

II 

oo 

o 

0.5329 

0.5326 

0.5302 

0.8004 

0.8061 

0.8032 

x=16.0 

0.2981 

0.2978 

0.2969 

0.3740 

0.3814 

0.3800 

x=17.0 

0.2840 

0.2838 

0.2829 

0.3691 

0.3649 

0.3635 

x=18.0 

0.2736 

0.2733 

0.2725 

0.3715 

0.3782 

0.3770 

x=19.0 

0.2668 

0.2667 

0.2660 

0.3807 

0.3752 

0.3740 

x=20.0 

0.2640 

0.2639 

0.2631 

0.3960 

0.4083 

0.4083 


Z 



A 

3.0 IN 


FIGURE 6 
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CONCLUSIONS 


The treatment of body forces of the centrifugal, gravitational, and thermal types in the 
implicit-differentiation formulation for the design sensitivity analysis of two-dimensional, 
axisymmetric, and three-dimensional problems is presented. The particular integral sensitivity 
expressions for the gravitational and centrifugal type body forces are developed. The thermoelastic 
sensitivity kernels are given for the thermal type body forces. A semi-analytical and a full-analytical 
approach for determining the sensitivity system matrices are used. A wide range of problems are 
solved for design sensitivities due to body forces and the results are validated through comparisons 
with exact analytical solutions. 
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